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ABSTRACT 

Reynolds  stress  closure  models  based  on  the  recursion  renormalization  group  theory  are  devel¬ 
oped  for  the  prediction  of  turbulent  separated  flews.  The  proposed  model  uses  a  finite  wavenum¬ 
ber  truncation  scheme  to  account  for  the  spectral  distribution  of  energy.  In  particular,  the  model 
incorporates  effects  of  both  local  and  nonlocal  interactions.  The  nonlocal  interactions  are  shown 
to  yield  a  contribution  identical  to  that  from  the  e-RNG,  while  the  local  interactions  introduce 
higher  order  dispersive  effects.  A  formal  analysis  of  the  model  is  presented  and  its  ability  to  accu¬ 
rately  predict  separated  flows  is  analyzed  from  a  combined  theoretical  and  computational  stand 
point.  Turbulent  flow  past  a  backward  facing  step  is  chosen  as  a  test  case  and  the  results  obtained 
based  on  detailed  computations  demonstrate  that  the  proposed  recursion  -  RNG  model  with  finite 
cut-off  wavenumber  can  yield  very  good  predictions  for  the  backstep  problem. 


t  This  research  was  supported  by  the  National  Aeronautics  and  Space  Administration  under  NASA  contract  No. 
NAS1-19480  while  the  authors  were  in  residence  at  the  Institute  for  Computer  Applications  in  Science  and  Engi¬ 
neering  (ICASE),  NASA  Langley  Research  Center,  Hampton,  VA  23681-0001. 


1.  INTRODUCTION 


Turbulent  flows  of  scientific  and  engineering  importance  are  characterized  by  a  broad  spectrum 
of  length  and  time  scales.  While  the  physical  aspects  of  turbulent  flows  are  best  described  by  the 
equations  of  motion,  limitations  in  computer  capacity  and  speed  preclude  their  direct  solution  for 
complex  flows  of  relevance  to  technical  applications.  The  current  practice  for  high  Reynolds  num¬ 
ber  flows  of  practical  interest  therefore  involves  some  type  of  modeling  for  Reynolds  stresses.  The 
commonly  used  turbulence  models  are  based  on  the  calculation  of  one-point  first  and  second  mo¬ 
ments  such  as  the  mean  velocity,  mean  pressure  and  turbulent  kinetic  energy.  Among  these,  the 
two-equation  turbulence  models  that  involve  the  use  of  transport  equations  for  the  turbulent  field 
parameters  that  involve  the  length  and  the  time  scales  are  probably  the  most  widely  used.  They  in¬ 
volve  the  simplest  level  of  Reynolds  stress  closure  that  do  not  depend  specifically  on  the  flow  ge¬ 
ometry.  (For  an  excellent  review  of  recent  trends  in  analytical  methods  for  Reynolds  stress  closure, 
the  reader  is  referred  to  Speziale.1) 

In  its  standard  form  the  two-equation  Reynolds  stress  turbulence  models  involve  the  turbulence 
kinetic  energy  and  dissipation  based  on  a  Boussinesq  type  approximation2  of  the  form 
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wherein  U  is  the  mean  velocity  based  on  Reynolds  average,  K  is  the  turbulence  kinetic  energy,  and 
V/'is  the  eddy  viscosity  which  is  isotropic.  Such  a  representation  of  turbulence  is  often  not  effective 
from  both  theoretical  as  well  as  phenomenological  point  of  view  and  the  shortcomings  associated 
with  it  are  fully  discussed  by  Speziale.1  To  over  come  some  of  these,  models  that  are  nonlinear 
(i.e.,  quadratic)  in  the  mean  strain  rate  were  proposed  in  the  form  of  a  constitutive  relation.3’4  Spe¬ 
ziale3  employed  tensor  and  dimensional  analysis,  together  with  invariance  constraints,  to  derive  his 
model.  Yoshizawa’s  model4  was  obtained  by  appealing  to  a  two-scale  direct  interaction  approxi¬ 
mation.  The  application  of  these  models  depend  on  the  empirical  evaluation  of  the  model  constants. 
Apart  from  the  specific  values  of  the  constants  in  these  models,  these  two  models  have  quite  similar 
structure  and  both  were  able  to  predict  the  anisotropy  in  the  Reynolds  stresses  in  a  noncircular  duct 
problem. 


The  present  study  addresses  the  need  for  a  more  effective  approach  in  the  development  of  two- 
equation  turbulence  models,  and  in  this  context  the  renormalization  group  ( RNG )  theory  based 
models  are  examined  for  further  development.  While  the  application  of  renormalization  group  the¬ 
ory  to  turbulence  has  attracted  much  attention,5'18  it  is  important  to  realize  that  these  calculations 
fall  into  two  distinct  categories:  (a)  e-RNG 5‘9,  pioneered  by  Forster  et.  a/.14,  and  (b)  recursion- 
/?A/G,10'13  pioneered  by  Rose.15  These  techniques  have  been  critiqued6,16'18  and  compared.12 


Here,  we  wish  to  point  out  that  in  the  e  -RNG,  a  small  parameter  €  is  introduced  into  the  exponent 
of  the  forcing  correlation  function  (with  the  forcing  function  being  introduced  into  the  momentum 
equation).  The  theory  is  then  developed  for  e  «  1,  and  all  constants  generated  are  evaluated  in 
this  limit  e  «  1.  However,  at  the  same  time,  all  exponents  that  are  e  -dependent  are  evaluated  at 
e  =  4.17  In  fact,  e  =  4  is  required  in  the  e  -RNG  to  recover  the  Kolmogorov  energy  spectrum  in 
the  inertial  range.17  This  not  only  plays  some  havoc  with  the  evaluation  of  constants,  it  also  leads 
to  another  problem:  /?lVG-induced  interactions  that  can  be  shown  to  be  irrelevant  in  the  limit 
e  «  1  cannot  now  be  shown  to  be  irrelevant  in  the  limit  e  — >  4.  Yet  in  e  -RNG  theories  (which 
require  s  — » 4)  these  higher  order  nonlinearities  are  assumed  to  be  unimportant.6,52,16  Moreover, 
e  -  RNG  theory  can  only  take  into  account  non-local  interactions.6’10’16  However,  the  procedure 
is  quite  amenable  and  Rubinstein  and  Barton7  have  derived  a  Reynolds  stress  model  using  €  -RNG 
methods  that  is  qualitatively  similar  to  that  of  Speziale3  and  Yoshizawa.4 

On  the  other  hand,  recursion-/?M7  does  not  rely  on  an  e  -expansion,  and  treats  explicitly  the 
cubic  nonlinearities  induced  into  the  renormalized  momentum  equation.  Moreover,  recursion-/?A?G 
can  handle  both  local  and  non-local  interactions.  Effects  such  as  the  cusp  behavior  in  the  transport 
coefficients  as  k  ->  kc  are  recovered  in  these  theories  (here,  kc  is  the  cutoff  wavenumber  separating 
the  large  scale  from  the  local  resolvable  scales).13,19'24  These  effects  are  the  consequences  of  local 
interactions  and  the  cubic  nonlinearities  introduced  by  the  RNG  procedure.  However,  one  of  the 
major  difficulties  to  the  application  of  recursion-/?M7  for  turbulent  flows  governed  by  the  Navier- 
Stokes  equation  was  that  the  eddy  viscosity  was  determined  as  a  fixed  point  of  a  very  complicated 
integro-dif fere  nee  equation.  This  drawback  has  now  been  removed  by  extending  the  theory  to  han¬ 
dle  the  iterative  removal  of  infinitesimal  wavenumber  bands.25  Now,  as  in  e-RNG,  the  eddy  vis¬ 
cosity  is  readily  determined  from  the  solution  of  a  relatively  simple  differential  equation.  Unlike 
e  -RNG,  though,  the  transport  coefficients  are  determined  over  the  whole  resolvable  scales  and  not 
just  in  the  wavenumber  limit  k  — » 0. 

In  the  present  work,  the  recursion  RNG  procedure  is  used  to  develope  a  Reynolds  stress  closure 
model  in  a  formal  manner.  It  is  shown  that  the  recursion-/?iVG  based  model  introduces  two  addi¬ 
tional  terms  arising  from  the  local  interactions  effects  that  are  of  higher  order  than  those  obtained 
using  the  e  -RNG  method  by  considering  only  the  long  wavelength,  non-local  interaction  limit.  The 
first  of  these  two  terms  as  well  as  those  from  the  conventional  turbulence  models  are  shown  to  be 
a  part  of  the  integrity  basis  used  in  the  representation  of  the  anisotropic  part  of  the  Reynolds  stress 
tensor.  °  The  second  term  which  arises  from  pressure- strain  coupling  is  quartic  in  strain  rate  and 
is  of  the  same  order  as  the  remaining  terms  of  the  integrity  basis. 

The  model  is  then  applied  for  turbulent  flow  past  a  backward-facing  step  which  has  played  a 
central  role  in  benchmarking  the  performance  of  turbulence  models  for  separated  flows.  During  the 


-2- 


past  decade  —  beginning  with  the  1980/81  Stanford  conference  on  complex  turbulent  flows29  — 
various  two-equation  turbulence  models  have  been  tested  and  compared  with  the  experimental  data 
of  Kim,  Kline  and  Johnston30  and  Eaton  and  Johnston31  for  the  backstep  problem.  Initial  results29 
indicated  that  the  standard  AT-e  model,  with  wall  functions,  underpredicted  the  reattachment  point 
by  a  substantial  amount  on  the  order  of  20-25%.  Several  independent  studies  have  been  subse¬ 
quently  published  using  alternative  forms  of  the  Af-e  model  wherein  a  variety  of  conflicting  results 
have  been  reported.  Considering  the  need  to  accurately  predict  separated  turbulent  flows  —  which 
can  have  a  wealth  of  important  scientific  and  engineering  applications  —  the  proposed  model  is 
applied  for  the  backstep  problem.  The  computations  based  on  a  sufficiently  resolved  finite-volume 
algorithm  show  that  the  proposed  model  based  on  the  recursive  application  of  the  renormalization 
group  theory  (developed  independently  without  any  ad  hoc  empiricism)  can  yield  a  prediction  for 
the  reattachment  point  that  is  within  a  few  percent  of  the  experimental  result.  The  physical  impli¬ 
cations  that  these  results  have  will  be  also  discussed  in  detail  in  the  sections  to  follow. 
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2.  FORMULATION  OF  THE  PHYSICAL  PROBLEM 


The  turbulent  motion  of  viscous,  incompressible  fluids  are  governed  by  the  Navier-Stokes 
equations  which  may  be  analyzed  using  classical  single  point  closure  based  on  Reynolds  decom¬ 
position  of  all  physical  variables.  The  resulting  averaged  equations  of  motion  are  of  the  form 


at/,  at/,.  _  _ a/>  a 2ui  _^ia 

dt  +  adxa  a*,. +  V°dxadxa  dxa 


at/a 

dxa 


=  o 


a) 

(2) 


where  t/;  is  the  mean  velocity,  P  is  the  mean  pressure,  Vq  is  the  kinematic  viscosity  of  the  fluid  and 
t:j  is  the  Reynolds  stress  tensor.  While  the  physical  aspects  of  turbulent  flows  are  best  described 
by  the  above  equations,  limitations  in  computer  capacity  and  speed  preclude  their  direct  solution 
for  complex  flows  of  engineering  importance.  The  current  practice  for  high  Reynolds  number 
flows  of  engineering  interest  therefore  involves  some  type  of  modeling  for  Reynolds  stresses.  The 
commonly  used  turbulence  models  are  based  on  the  calculation  of  one-point  first  and  second  mo¬ 
ments  such  as  the  mean  velocity,  mean  pressure  and  turbulent  kinetic  energy. 

In  the  present  work  we  consider  the  development  of  Reynolds  stress  model  by  the  recursion 
RNG  formulation  and  its  application.  In  this  context  it  is  convenient  to  express  the  equations  of  mo¬ 
tion  using  the  Founer  representation, 

(^  +  vo*2)w/(k>f)  =  jd3P  «a(P.0Kp(k-p,r)  (3) 

*aMk’  0=0  (4) 

It  is  important  to  note  here  that  no  random  forcing  is  introduced  unlike  in  the  €  -RNG  theories 
where  it  plays  a  critical  role  in  introducing  the  small  parameter  € .  The  nonlinear  coupling  coeffi¬ 
cient 


Mia p(« 

with  Dap  being  the  projection  operator  defined  by 


(5) 


D„|S  (*>  =  " 


Vis 


an  k2 


(6) 


A  scale  factor  h  is  now  introduced  to  partition  the  wavenumber  space  to  N  segments  such  that 
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kc  =  kN  =  hNkQ  , 

*jv_i=  hN~lk0  ,  ... ,  kx  =  hkQ  (7) 

k0  ~  0(kd) 

where,  kc  is  the  cut-off  wavenumber  that  separates  the  large  scale  from  the  small  scale,  and  kd  is  of 
the  order  of  Kolmogorov  dissipation  wavenumber.  The  velocity  field  is  then  decomposed  in  the 
Fourier  space  such  that  ui  =  uf  0  (£  -  kh)  +  uf  0  {kh  -  k)  (wherein  0  is  the  Heavyside  unit  step 
function,  k /,  j  is  the  local  cut-off  wavenumber,  the  superscript  <  represents  large  resolvable 
scale  quantities  while  >  represents  small  unresolved  scales  with  respect  to  the  local  cut-off  wave- 
number).  The  evolution  of  these  components  can  be  directly  obtained  from  (3)  as: 

(jt  +  V2)  uf  (k,  t)  =  A/ia(J  (k)  jd3p  [u*( p,  t)  (k  -  p,  t) 

+  2u*  (p,  r)  ujj  (k  -  p,  t)  +  ul  (p,  t)  up  (k  -  p,  t )] 

(5r +  vo*2) (k>  =  jd3P  (P»  0  «p  (k  ~  P-  0 

+  2u>(p,  r)up  (k-p,  t)  +u>  (p,  r)  uj  (k-p,  t)] 

It  should  be  noted  that  in  e-RNG  approach,  one  is  forced  into  taking  the  large-scale  infrared 
limit  k  -*  0.  In  essence,  this  forces  a  spectral  gap  between  the  resolvable  part  of  the  flow  field  and 
the  small  unresolved  scales.  If  this  spectral  gap  were  somehow  present  initially,  it  would  be  quickly 
populated  in  just  a  few  eddy  turn  over  times.  Thus,  retaining  only  the  distant  interactions  may  not 
be  appropriate.  In  fact,  it  has  been  shown13  that  the  energy  transfer  function  that  corresponds  to 
local  interactions  accounts  for  most  of  the  energy  flow  out  of  the  resolvable  scales.  It  thus  seems 
important  to  retain  both  local  and  nonlocal  interactions  in  the  modeling  of  the  Reynolds  stress  and 
this  can  be  readily  achieved  by  recursion-/WVG.  In  particular,  it  is  apparent  that  the  Reynolds  stress, 
x.j  =  if?  +  xff ,  has  two  components.  The  if*  -  part  arises  from  the  infrared  limit  of  k  — >  0  and 
is  due  to  the  uf  -  uf  distant  interaction  limit  while  the  xf?  -  part  arises  from  the  0  <  k  <  kc  spec¬ 
trum  and  is  due  to  the  uf  -  uf  local  interaction  limit.  Thus  in  the  e  -RNG  model,  the  Reynolds 
stress  =  xf?  and  is  obtained  purely  from  the  uf  -  uf  interaction  in  the  small  unresolved  scale 
momentum  equation  (8).  We  now  consider  the  contribution  to  the  Reynolds  stress  that  arises  from 
the  local  interaction: 


XV  ~  ~jd3P[u?  (k  ~  P)  u>j  (P>  +  uj  <k  ~  P)  (P)  1  (10) 

Herein  u<(k-p)  corresponds  to  the  Fourier-transformed  velocity  field  in  the  resolvable  large  scales, 
Ik-pl  <  kv  while  u>(p)  corresponds  to  the  small  scale  field  with  Ip!  >  and  kl  is  the  wavenumber 
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which  separates  the  resolvable  from  the  small  unresolved  scales. 

Now  there  are  certain  constraints  that  t?5  must  satisfy;  in  particular, 

(a)  — >  0  as  the  turbulent  kinetic  energy  K  -»  0. 

(b)  is  Galilean  invariant. 

Constraint  (a)  is  clear:  as  the  turbulent  kinetic  energy  K  — >  0,  the  small  scale  velocity  field  iC  -4  0, 
so  that  — >  0.  Constraint  (b)  arises  since  it  has  been  shown  in  detail  that  recursion-/?M7,  even 

with  its  cubic  nonlinearity  interaction,  is  a  Galilean  invariant  theory.25  Now,  if  one  follows  stan¬ 
dard  recursion-ZWG  procedures,10'13  the  relevant  part  of  the  small  scale  velocity  field  that  contrib¬ 
utes  to  TtJc  is  (on  the  removal  of  the  first  subgrid  shell,  c.f.,  equation  (9)). 

A/jcjQ  (p)  m 

“i  (p)  = - r~rJ'  u«(J)up  0  —  p)  +  •••  (ii) 

V0 P 

where, . . .  refers  to  terms  that  will  not  contribute  to  x** .  At  the  removal  of  the  first  small  scale  shell 
Vq  is  the  molecular  viscosity.  At  subsequent  steps  in  the  removal  of  the  subgrid  scales,  v0  will  be 
replaced  by  the  renormalized  eddy  viscosity. 

In  the  substitution  of  the  subgrid  velocity  field  from  (11),  into  the  Reynolds  stress  given 
by  (10),  one  must  ensure  that  the  constraints,  (a)  and  (b)  above  are  satisfied.  It  can  be  shown  that 
these  constraints  are  not  satisfied  by  a  simple  substitution  of  (1 1)  imo  (10).  Hence  we  introduce  a 
factor  (p  —  k)a(p  —  s)a/k*  so  as  to  recover  these  properties.  This  factor  is  appropriate  since, 
from  (11),  kc  =  0(dK3/2)  — » oo  as  K  — >  0.  Hence,  in  recursion-RNG,  after  the  elimination  of  the 
first  spectral  band  (i.e.,  after  the  first  iteration)  (10)  becomes 


oo .  =  -4 


(Pa  ^a)  ( Pa  ^a) 


v  n* 


Mjys  (p)  U *  (s)  Mg  (p  -  s )  uf  ( k-p )  +  i 


(12) 


where  p  is  in  the  subgrid  region,  while  the  arguments  of  all  the  velocity  fields  are  in  the  resolvable 
scales.  The  second  term  in  (12)  is  obtained  by  i  j  interchange  in  the  first  term.  At  the  next  ap¬ 
plication  of  the  above  procedure  (i.e.,  after  the  second  iteration) 


T5<(k)i2  =  t><(k)|rp3p[W<(k-P)M;  (P)  +u<  (k  —  p)  n(>  (p)j  (13) 


where,  k  <  ^2 •  For  k2<k<  fcj,  the  small  scales  are  obtained  by  a  similar  procedure  by  utilizing  the 
equations  of  motion  for  small  scales.10,15  The  resulting  expression  for  the  small  scales  may  then 
be  written  as: 


(P) 


— — — — 1  fiV 

V,  (p)p  J 


«a  (j)«J  (j-P)+- 


(14) 
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substituting  from  (12)  and  (14)  in  (13), 


tr  (k)  = 


V  f  C  .3  j3  (Pa~~  <  <  < 

X  - Mjy&(p)uy  (s)ub  ( P-S)ui  (k-p)  +  t «-»/ 


'ch=  iL 


v*  (p)pl 


Proceeding  similarly  to  the  elimination  of  the  N-th  spectral  band,  one  finds 
rX 1  V*  T  fw3  v3  (PoL~^a)  (P  a~so) 

^  W  -  --2  £ 


veA=  1  L 


^Mp 


(15) 


Mjyb(p)uy  (j)wg  (p-s)uf  (k-p)  +  /  <->y 


(16) 


where,  £c  is  the  cut-off  wavenumber  that  separates  the  large  scales  from  the  small  scales.  In  (16) 
the  summation  is  over  all  the  spectral  bands  eliminated  and  V/,(p)  is  the  renormalized  eddy  viscos¬ 
ity  . 10-13  We  now  proceed  to  the  differential  limit  of  infinitesimal  spectral  shells,  so  that15 


1 

v(*c) 


(17) 


for  a  Kolmogorov  /:  5j^ 


inertial  range  energy  spectrum.  Since  the  quadratic  coupling  coefficient 


-  s(V*  + Vf"2^) 

in  the  differential  limit  for  infinitesimal  spectral  shells  of  the  small  scales,  (16)  becomes 


(18) 


(k)  =  - 


1 


2  tv  (*,)  k 


10/3 


f  j3  j3  0>«  (^o  ^ 

jd3pd3s - 271 - (Pfiuj  (/>-*)«,■  (k-p) 


+  p$u<(s)uf(p-s)u*(k-p)  -2  (s)  uy  (p  -  s)  uf  (k- p)  }  +i<r+j 


(19) 

The  inverse  transform  of  (19)  yields  the  contribution  to  Reynolds  stress  due  to  local  interaction  in 
physical  space  and  is  shown  below. 


Tj*  (x) 


(2tc)M*c)*‘ 


10/3 


dU^x)  f  ( 


dX'r 


+  g . ..  a  f^coa2^*') 

2  dx'i  dx'adx'y  j 


dx' 


}  + '  <~>j 


(20) 


In  the  above,  the  resolvable  velocity  field  in  the  physical  space  is  identified  with  the  mean  velocity 
£/,•>  and  the  nonlocal  kernels  -  x’)  and  g2(x  -  x')  are  defined  by: 
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„  ,_N  _  [exP  (ip '  r)  ^ 

r  J  p2/3  ^ P 


and 


rexp  (ip  ■  r)  , 

«2<r>  °  J  -'i/r 


<21. 


The  structure  of  these  nonlocal  kernel?  > ■>  dependent  on  the  cut-off  wavenumber,  kc.  Now, 

2kr  2 


gl(r)  -  ~  j  /?1/3sin  (pr)dp  =  ^**/3JV  1/3sin  (kcp'r)  dp' 


An 


(22) 


l 

where  the  integration  limits  are  obtained  from  the  partial  average1^  over  the  nearest  unresolved 
band.  Similarly, 


g2(r)  =  ~r *c2/3jV”5/3sin  {kcp'r)dp'.  (23) 

t 

These  nonlocal  kernels  reduce  to  local  expressions  for  kc  »  1 .  Noting  that  sin  ( kLJ. :)  /x  — >  rc5  (x) 
as  kc  oo, 


gt(r)  /4jt  =  3.58^/36(r) 

g2(r)/4tt*  I.74Jfc"2/3d(r) .  (24) 

Following  Yakhot  and  Orszag,5  this  cutoff  wavenumber  kc  can  be  related  to  the  isotropic  part  of 
the  turbulent  kinetic  energy  K  as: 


K  = 


80 


3  e2/3 

”  2C**2/3 
*c 


(25) 


on  using  the  Kolmogorov  internal  range  spectrum  for  E(p),  and  where  Cg  is  the  Kolmogorov  con¬ 
stant.  Thus,  the  cut-off  wavenumber 


e  3  3/2 

* C  ~  2 

Hence  under  this  approximation  Tj><  (x)  in  (20)  reduces  to  a  simple  algebraic  form 


(25) 


=^1^3 


dU.dU  ,dUn 
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+  C 


K 

R27 


dxa  dxy  dxa 
where  use  has  been  made  of  Kraichnan’s'  ^  result 

v{kc)k*/3  =  0.19CJe1/3 


dUl  d 

dxt1  dx 


fdUp  d2Uy 

fit.,  f).t,.fi.tn 


+  i  *-»  / !  (27) 


(28) 
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However,  for  most  flows  of  interest,  kc  varies  considerably  throughout  the  flow  domain.  Thus,  if 
one  wishes  to  use  the  local  approximation  for  the  kernels,  gj  and  gi  as  well  as  retain  kc~  0(1) 
effects  so  as  to  obtain  an  algebraic  form  for  x^  (x) ,  then  the  coefficients  and  Crj  will  be 
function  of  the  flow  quantities.  In  particular,  Crj  and  C may  be  expressed  as: 

CRl  ==  Cri(e,K,^0  and  ^R2  ~  ^R2  C)  (29) 


where 


the  strain  rate,  r\  s  (SapSap) 1/2 ,  with  Sap  =  - 


3*0  dxa  , 


and 


j  fdUa  5(/p\ 

the  rotation  rate,  C  s  (WapWap)  U2>  with  Wafi  -  \  (^3^  "  dx~  J' 


Rubinstein  and  Banon7  have  derived  a  Reynolds  stress  model  using  the  6  -RNG  method  (which 
corresponds  to  the  infrared  limit  of  k  — >  0)  of  the  following  form: 


2  c  dui 

%K  5.  +  vt  — 

3  v 


dUA 


■WtWj}* 

dxadxa) 


+  Cx2 


(dU^Ua 
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dUjdUa\* 
dxa  dx-  j 


+  cx3 


*  -J 

[dXi  dXj  J 

(30) 


The  constants  Ctl  =  0  .034,  Cx 2  =  0. 104  and  Cx 3  =  -0.014,  and  (...)*  denotes  the  deviatoric  part  of 
the  expression  within  the  parenthesis.  The  first  two  terms  correspond  to  the  linear  model  and 
vr  =  is  the  isotropic  eddy-viscosity,  where  e  is  the  turbulence  dissipation  and  =  0.09 

based  on  empirical  data  from  equilibrium  boundary  layer  flows.  The  above  model  is  quadratic  in 
mean  strain  rate,  includes  the  effect  of  convection  and  diffusion  and  is  qualitatively  similar  to  other 
second  order  models  3,4 


Combining  (27)  -  (30)  a  formal  expression  for  Reynolds  stress  which  includes  both  the  local 
and  nonlocal  interactions  may  be  expressed  in  the  following  form; 
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(31) 


It  is  of  some  interest  to  recast  the  above  expression  for  the  Reynolds  stress  obtained  from  the 

recursion  RNG  model  into  an  integrity  basis  representation  which  are  commonly  employed  to  rep- 

2 

resent  the  anisotropic  part  of  the  Reynolds  stress  tensor,  btj  =  ( t  - .  -  ^Kb . .)  /2 K  in  the  following 
form-invariant  manner. 


T“ 


where,  T-p  is  the  integrity  basis  for  functions  of  symmetric  and  antisymmetric  tensor  and  the  co- 
efficients  G  are  scalar  functions  of  the  irreducible  invariants  of  the  strain  rate  tensor,  5,,  and  the 
rotation  rate  tensor,  Wy 26-28 


T(i)  _  c 
Jij  ~  *ij 

=  C  W  -W  S 
1  ij  0 iawaj 

-tO)  _  c  c  _  1  C  X 
1  ij  JiaJcy  g^aa0!; 

TiP  =  V,- 

rf  =  wias^sfj-  s,asa6  wfj 

PP  =  Vopsw + -  \  (s^iw  8aP 

rf  =  +  VanV-u  -  3  (S.pVteW'c.)  6„ 

4I0>  =  (33) 


Recently,  Gatski  and  Speziale28  have  applied  this  integrity  basis  to  determine  an  explicit  alge¬ 
braic  stress  model  for  three  dimensional  turbulent  flows  based  on  a  systematic  derivation  from  a 
hierarchy  of  second-order  closure  models.  In  their  applications,  they  have  restricted  themselves  to 
a  model28  which  involves  just  the  tensors  T-^  ,  t\P  ,  and  T-p  .  It  can  readily  be  shown  that  the 
k  =  0  part  of  the  Reynolds  stress,  x~> ,  specified  in  (30),  involves  the  tensors  Tp ,  T-p  ,  T-p  and 
T-p  from  the  integrity  basis  given  by  (33).  The  z??  -part  of  the  Reynolds  stress  is  common  to  both 
6  -RNG  andrecursion-RNG.  However,  unlike  e  -RNG,  there  is  now  a  finite-^  spectral  contribution, 
XyK  to  the  Reynolds  stress  in  recursion-/?ZVG.  It  can  be  shown  that  the  first  term  in  (29)  involves 
the  tensors  T^p  and  T\p  of  the  integrity  basis  (33).  The  second  term  in  (29),  however,  is  not 
readily  expressed  in  this  second  order  basis  because  of  the  intrinsic  fourth  rank  tensors  involved. 


-  10- 


However,  as  can  be  seen,  this  term  is  quartic  in  strain  rate  and  arises  from  pressure-strain  coupling, 

f7)  (8) 

and  could  be  thought  of  as  having  contributions  similar  to  those  involving  the  tensors  7-  ,  7--  , 

and  Tjp . 

The  above  expression  for  Reynolds  stress  (31)  are  to  be  used  along  with  the  equations  of  mo¬ 
tion  (1)  by  specifying  turbulent  kinetic  energy  and  dissipation.  In  two-equation  turbulence  models, 
this  closure  is  achieved  through  the  development  of  transport  equations  for  the  turbu:  nt  kinetic 
energy  and  dissipation  —  quantities  that  are  directly  related  to  the  length  and  time  scales  —  of  the 
following  general  form: 


3  K  3  K  ^  3 

-3—  +  (A'5 —  =  rP~  £  +  -5 — 
3 1  Jdx;  dx; 


0  aKJdxi 


3e  3e  „  e  _  „  e2  3 
_  +  f/  _  Cel-?P-Ce2-  +  ^- 


V°+«j3xt 


where,  v  =  Vq  +  Vj  is  the  total  viscosity,  IP  =  (dU/dXj)  is  the  turbulence  production,  e  is  the 
scalar  turbulent  dissipation  rate,  and  2v0S,y  (3w /3x()  (3u-/3xy)  is  the  turbulent  strain  rate 
correlation  term.  The  quantities  CE;,  C&,  ae  are  dimensionless  and  taken  to  be  1.44,  1.92,  1.0 
and  1.3,  respectively,  consistent  with  the  standard  form  of  the  two-equation  K-e  model  (based  on 
empirical  data  obtained  from  equilibrium  boundary  layer  flows).4  It  should  be  noted  here  that  these 
are  calculated  explicitly  to  be  1.42,  1.68,  0.719  and  0.719,  respectively,  in  the  e  -ftVG-based  for¬ 
mulations.5,6  In  addition,  it  is  customary  to  neglect  the  strain  rate  correlation  term  3(,in  the  standard 
K-e  model  as  well  as  in  some  RNG-based  models7  since  3£must  be  modeled  in  order  to  achieve 
closure.  Special  forms  of  approximations  have  been  recently  employed  in  the  case  of  some  of  the 
e  -/?M7-based  models.8,9  To  avoid  ambiguity  (and  to  avoid  introducing  ad  hoc  assumptions  re¬ 
quired  to  model  the  contribution  from  the  turbulent  strain  rate  correlation  term  is  neglected  in 
the  present  work. 
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3.  TURBULENT  FLOW  PAST  A  BACKWARD-FACING  STEP  —  A  CASE  STUDY 

The  problem  to  be  considered  is  the  fully-developed  turbulent  flow  of  an  incompressible  vis¬ 
cous  fluid  past  a  backward-facing  step  (a  schematic  is  provided  in  figure  1).  Calculations  will  be 
conducted  for  an  expansion  ratio  (step  height:  outlet  channel  height)  E  of  1:3  and  the  Reynolds 
number  Re  =  132,000  based  on  the  inlet  centerline  mean  velocity  and  outlet  channel  height  (which 
corresponds  to  that  of  Kim  etal?Q  and  Eaton  and  Johnston31).  The  mean  turbulence  equations  (cf., 
§2)  are  solved  subject  to  the  following  boundary  conditions:33 

(a)  inlet  profiles  for  U,  K  and  e  are  specified  five  step  heights  upstream  of  the  step  comer  (U  is 
taken  from  the  experimental  data30'31  and  the  corresponding  profiles  for  K  and  £  are  computed 
from  the  model  formulated  for  channel  flow), 

(b)  The  law  of  the  wall  is  used  at  the  upper  and  the  lower  walls,  and 

(c)  Conservative  extrapolated  outflow  conditions  are  applied  thirty  step  heights  downstream  of  the 
step  comer;  these  conditions  involve  the  following:  i)  the  V'-component  of  the  velocity  for  the 
cells  at  the  outflow  boundary  are  obtained  by  extrapolation;  ii)  the  (/-component  of  the  velocity 
is  then  computed  by  the  application  of  a  mass  balance;  and  iii)  the  scalar  quantities  such  as  pres¬ 
sure,  turbulent  kinetic  energy  and  turbulent  dissipation  are  all  obtained  by  extrapolation.  It  was 
found  that  a  downstream  channel  length  of  about  thirty  step  heights  was  needed  to  ensure  that 
the  local  error  for  all  the  quantities  was  of  the  same  order  as  the  interior  values. 

A  finite  volume  method  which  relies  on  solving  the  discretized  equations  by  a  line  relaxation 
method  with  the  repeated  application  of  the  tridiagonal  matrix  solution  algorithm  is  modified  and 
applied  for  the  present  case  to  obtain  the  steady  state  solution.33,34  The  computed  solution  was  as¬ 
sumed  to  have  converged  to  its  steady  state  when  the  root  mean  square  of  the  average  difference 
between  successive  iterations  was  less  than  10'4  for  the  mass  source.34  Approximately  2000  itera¬ 
tions  were  needed  for  the  convergence  of  the  standard  K-z  model;  this  corresponds  to  approximate¬ 
ly  20  minutes  of  CPU  time  in  a  partially  vectorized  mode  on  the  Cray-YMP  supercomputer  using 
64-bit  precision.  The  recursive  fttVG-based  K-z  model  requires  approximately  33%  more  CPU  time 
due  to  the  fact  that  the  additional  terms  in  the  RNG- based  K-z  model  have  to  be  evaluated  during 
each  iteration.  These  correction  terms  are  dispersive  —  an  additional  feature  that  slows  conver¬ 
gence. 

The  issue  of  resolution  is  crucial  for  the  backstep  problem  and  calculations  indicate  that  a 
200x  1 00  mesh  yields  a  fully  grid  independent  solution.33  All  of  the  computations  conducted  in  this 
study  were  performed  using  this  200x100  nonuniform  mesh.  As  indicated  earlier,  the  inlet  condi¬ 
tions  were  specified  5  step  heights  upstream  of  the  step  comer  and  the  outlet  boundary-  conditions 
were  specified  30  step  heights  downstream  of  the  step  corner.  It  is  crucial  that  a  sufficient  distance 
downstream  of  the  reattachment  point  be  allowed  before  the  outflow  conditions  are  imposed.  Many 
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earlier  computations  of  the  backstep  problem  were  in  significant  error  due  to  the  imposition  of  ful¬ 
ly-developed  outflow  conditions  too  close  to  the  reattachment  point.  Furthermore,  it  is  crucial  that 
a  fine  mesh  be  used  near  the  step  comer  for  computational  accuracy.  It  should  be  also  noted  that 
the  law  of  the  wall  does  not  formally  apply  to  separated  turbulent  boundary  layers.  However,  since 
the  separation  point  is  fixed  at  the  corner  of  the  backstep  —  and  the  flowfield  is  solved  iteratively 
so  that  the  friction  velocity  ux  can  be  updated  until  it  converges  —  major  errors  do  not  appear  to 
result  from  its  use. 

First,  results  will  be  presented  for  the  standard  K-e  model.  For  this  case  —  as  well  as  the  other 
results  to  follow  —  computed  results  for  the  mean  velocity  streamlines,  the  streamwise  mean  ve¬ 
locity  profiles,  the  streamwise  turbulence  intensity  profiles  and  the  turbulence  shear  stress  profiles 
are  compared  with  the  Kim  et  al?°  experimental  data  as  updated  by  Eaton  and  Johnston.31  In  fig¬ 
ure  2(a)  the  computed  streamlines  are  shown,  indicating  reattachment  at  Xr'H  =  6.12  —  a  result 
which  is  approximately  15%  underprediction  of  the  experimental  reattachment  point  of  X^H  =  7.1. 
In  figure  2(b),  the  streamwise  mean  velocity  profiles  predicted  by  the  standard  K-e  model  are  com¬ 
pared  with  the  experimental  data.  Except  in  the  vicinity  of  the  reattachment  point,  the  comparisons 
are  fairly  good.  More  serious  discrepancies  between  the  model  predictions  and  the  experimental 
data  occur  in  the  initial  part  of  the  recovery  zone  for  the  streamwise  turbulence  intensity  profiles 
as  shown  in  figure  3(a).  However,  the  model  predictions  for  the  turbulence  shear  stress  profiles  are 
reasonably  good  as  can  be  seen  from  figure  3(b). 

Next,  we  consider  the  computations  based  on  the  recursion -RNG  model.  For  computational  ef¬ 
ficiency,  the  nonlocal  kernels  gj  and  g2  are  approximated  by  one-dimensional  delta  functions.  This 
is  appropriate  if  kc  »  1 .  In  Fig.  4,  the  variation  of  kc  (normalized  by  the  step  height  and  expressed 
as  kcH)  is  shown  at  various  locations  downstream  of  the  step.  As  can  be  seen,  the  magnitude  of  kc 
is  at  least  an  order  of  magnitude  larger  than  the  grid  size  employed  {k^  >  1  for  the  200  x  100  mesh 
used).  Thus  the  local  representation  of  g\  and  gj  is  justified.  The  explicit  functional  form  of  the 
coefficients  Cg]  and  Cg2>  (29),  is  much  more  difficult  to  evaluate.  As  a  first  attempt  at  applying 
this  recursion-/?NG  model  we  make  a  lowest  order  approximation  that  these  coefficients  are  con¬ 
stants.  Clearly,  especially  in  light  of  the  great  success  in  enforcing  a  functional  ansatz  on  these 
types  of  coefficients,9  this  will  underestimate  the  robustness  and  accuracy  of  the  recursion-/?NG 
formulation.  Nevertheless,  it  is  deemed  more  appropriate  here  to  apply  the  theory  without  adding 
on  these  extraneous  functional  ansatzes.  It  is  also  evident  that  in  applying  the  one-dimensional  lo¬ 
cal  form  of  g,  (|x  -  x  j )  and  g2  (|x  —  x  j) ,  the  angular  dependences  of  the  second  term  in  (22)  is 
much  stronger  than  that  of  the  first  term  in  (22).  This  is  because  of  the  higher  order  derivatives 
present.  These  angular  integrations  are  also  dependent  on  spatial  position  due  to  the  changing  flow 
characteristics.  Again,  for  simplicity,  we  approximate  these  angular  integrations  by  constants,  and 
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thus  put  the  model  to  an  even  stricter  test.  In  particular,  for  the  present  analysis,  CR]  and  CR2  are 
taken  to  be  0.025  and  0.0342  x  10~3,  respectively. 

Now,  it  will  be  demonstrated  that  the  use  of  the  proposed  recursive  RNG- based  model  can  yield 
a  more  significant  improvement  in  the  results.  The  computed  streamlines  for  the  flow  field  shown 
in  figure  5(a)  have  a  mean  reattachment  point,  XrlH  »  6.72,  a  result  which  is  about  5%  lower  than 
the  experimental  value.  The  corresponding  mean  velocity  profiles  shown  in  figure  5(b)  indicate 
very  good  agreement  with  the  experimental  results.  It  should  be  noted  that  computations  performed 
based  on  the  e  -RNG- based  model  for  the  same  flow  conditions  and  model  constants  yield  a  mean 
reattachment  point,  X^H  ~  6.42  (not  shown  herein).  The  difference  in  the  size  of  the  separated  flow 
region  can  be  clearly  attributed  to  the  contributions  from  the  Reynolds  stress  terms  representing 
the  local  interactions  effects  in  the  recursion  RNG  model. 

Furthermore  the  overall  agreement  between  the  turbulence  intensity  and  the  shear  stresses  with 
the  experimental  data  shown  in  figures  6(a)  and  6(b)  are  also  good.  The  most  notable  difference 
between  the  predictions  of  the  RNG- based  models  and  the  standard  K-e  model  lies  in  the  stream- 
wise  turbulence  intensity  (figure  6(a)).  The  slight  trough  shaped  variation  predicted  in  this  region 
is  consistent  with  more  recent  independent  experiments.35 

In  addition,  to  illustrate  the  differences  associated  with  the  modeling  of  the  Reynolds  stresses, 
the  variation  of  the  turbulent  eddy-viscosity  normalized  with  respect  to  its  molecular  counterpart, 
v7-/v0,  is  shown  at  several  locations  downstream  of  the  backward-facing  step.  As  can  be  seen,  in 
the  recirculation  region  the  eddy  viscosity  predicted  by  the  standard  K-e  model  is  generally  larger 
than  that  due  to  the  recursion  RNG  leading  to  substantial  reduction  in  the  size  of  the  separated  flow 
region  (cf.,  figures  2(a)  and  5(a)). 

The  wall  pressure  coefficient  is  an  important  parameter  for  engineering  applications.  In  figures 
8(a)-(b),  the  pressure  coefficients  Cp  (=  2[p  -pr]/Ur2,  where  pr  and  Ur  are  the  reference  pressure 
and  velocity  which  are  taken  at  the  centerline  of  the  inlet)  obtained  from  the  standard  and  the  con¬ 
ventional  and  the  RNG-based  K-e  models  at  the  top  and  bottom  walls  are  compared  with  the  exper¬ 
imental  data  of  Eaton  and  Johnston.31  As  can  be  seen,  both  the  standard  and  the  RNG-based  K-e 
models  perform  comparably  well  in  reproducing  the  experimental  trends.  The  skin  friction  coeffi¬ 
cients  Cf~  2ux2/Ur2  obtained  from  the  standard  and  recurs  ion -RNG  K-e  models  are  compared 
with  the  scaled  experimental  data  of  Driver  and  Seegmiller36  for  the  bottom  wall  in  figure  8(c). 
Here,  we  make  use  of  the  fact  that  the  ratio  CftCf0 OI  when  taken  as  a  function  of  the  normalized 
distance  (X-Xr)/Xr,  is  independent  of  the  expansion  ratio  (given  that  Cfx  is  the  fully-developed 
skin  friction  coefficient  andXr  is  the  reattachment  point).  As  can  be  seen,  the  recursive  RNG- based 
K-e  model  performs  better,  however,  both  models  are  probably  within  the  uncertainty  of  the  exper¬ 
imental  data. 
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4.  CONCLUSIONS 


The  recursion  renormalization  group  ( RNG )  theory  is  utilized  to  develop  Reynolds  stress  clo¬ 
sure  models  for  the  prediction  of  turbulent  separated  flows.  Unlike  the  small  parameter  e  -RNG 
model  that  only  consider  the  interaction  of  the  long  wavelength  modes,  the  recursion-R/VG  model 
includes  both  the  local  and  nonlocal  interaction  of  all  the  relevant  resolvable  scales.  A  formal  de¬ 
velopment  of  the  Reynolds  stress  model  is  presented  and  the  resulting  higher  order  terms  by  con¬ 
sideration  of  the  local  interaction  effects  are  shown  to  be  qualitatively  similar  to  those  recently  ob- 
tained  by  considering  a  hierarchy  of  second-order  closure  models. 

The  ability  of  the  proposed  model  to  accurately  predict  separated  flows  is  analyzed  from  a  com¬ 
bined  theoretical  and  computational  standpoint  by  considering  turbulent  flow  past  a  backward  fac¬ 
ing  step  as  a  test  case.  The  results  obtained  based  on  detailed  computations  demonstrate  that  the 
proposed  recursion-&VG  model  can  yield  very  good  predictions  for  the  turbulent  flow  of  an  in¬ 
compressible  viscous  fluid  over  a  backward- facing  step.  It  should  be  remembered  that  the  deficien¬ 
cies  of  two-equation  models  are  well  established,  particularly  in  turbulent  flows  with  body  forces 
or  Reynolds  stress  relaxation  effects.1  Consequently,  the  findings  of  this  study  should  not  be  inter¬ 
preted  as  an  unequivocal  endorsement  of  two-equation  RNG  models.  Nonetheless,  this  study 
shows  that  properly  calibrated  two-equation  turbulence  models,  which  account  for  the  anisotropy 
of  the  turbulent  stresses,  can  be  effective  for  the  prediction  of  turbulent  separated  flows. 
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Figure  2  Computed  mean  flow  field  for  the  standard  K-e  model 
[E  =  1:3;  Re  =  132,000;  200x100  mesh) 
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Figure  3  Computed  turbulence  stresses  for  the  standard  K-z  model  [E  =  1 :3; 
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Figure  4  Computed  profiles  of  the  cut-off  wavenumber  for  the  recursion 
RNG  K-e  model  [E  =  1:3;  Re  =  132,000;  200x100  mesh;  —  computations; 
o  experiments  of  Kim  et  al,  1980;  Eaton  &  Johr  ston,  1981] 
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(b)  Dimensionless  mean  velocity  profile 
t - Computations; 

o  Experiments  of  Kim  etal,  1980;  Eaton  &  Johnston,  1981) 


Figure  5  Computed  mean  flowfield  for  the  recursion  RNG  /C-e  model 
[E  =  1 :3;  Re  =  132,000;  200x100  mesh] 
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Figure  6  Computed  turbulence  stresses  for  the  recursion  RNG  K-e  model 
[E  =  1:3;  Re  =  132,000;  200x100  mesh;  —  computations;  o  experiments  of 
Kim  ef  a/,  1980;  Eaton  &  Johnston,  1981] 
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Figure  7  Computed  profiles  of  eddy-viscosity  [E  =  1:3;  Re  =  132,000; 
200x100  mesh; - computations  for  the  standard  AC-e  model; - com¬ 

putations  for  the  recursion  RNG  K-e  model] 
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(c)  Skin  friction  coefficient  along  the  bottom  wall 

Figure  8  Comparison  of  the  predicted  wall  distribution  with  experiments 

I . computations  for  the  standard  K-z  model; - computations  for  the  re 

curslon  RNG  K-z  model;  o  experiments;  E  =  1:3;  Re  =  132,000] 
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